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ABSTRACT 

We live in a time where climate models predict future increases in environmen- 
tal variability and biological invasions are becoming increasingly frequent. A key 
to developing effective responses to biological invasions in increasingly variable en- 
vironments will be estimates of their rates of spatial spread and the associated un- 
certainty of these estimates. Using stochastic, stage-structured, integro-difference 
equation models, we show analytically that invasion speeds are asymptotically nor- 
mally distributed with a variance that decreases in time. We apply our methods to 
a simple juvenile-adult model with stochastic variation in reproduction and an il- 
lustrative example with published data for the perennial herb, Calathea ovandensis. 
These examples buttressed by additional analysis reveal that increased variability 
in vital rates simultaneously slow down invasions yet generate greater uncertainty 
about rates of spatial spread. Moreover, while temporal autocorrelations in vital 
rates inflate variability in invasion speeds, the effect of these autocorrelations on the 
average invasion speed can be positive or negative depending on life history traits 
and how well vital rates "remember" the past. 
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Introduction 



Anthropogenic forces are changing the temporal distribution of environmental fluctuations 
and accelerating the rate at which species are being introduced into non-native habitat. Gen- 
eral circulation models predict increased variability in temperatures and precipitation and these 



20 changes are likely to increase variability in vital rates of many species (Easterling et al. 2000 



Tebaldi et al. 2006 Boyce et al. 2006). Alternatively, human activities such as agriculture, 



Recreation, and transportation are spreading species beyond their natural dispersal barriers (El- 
1958; |Kolar and Lodge 2001). While most of these accidental or intentional introductions 



ton 



24 fail, the successful invaders can have devastating impacts on human health and native biodiver- 



25 sity (Kolar and Lodge 2001). To manage these impacts, it is essential to understand the rate 



of range expansion-the invasion speed-of these invaders. Here, we bring together the theory of 
stochastic demography and invasion speeds to provide a general framework to estimate invasion 
speeds for structured populations in a variable environment. 

Stochastic demography is concerned with understanding population growth and distribution 



30 
31 



when vital rates vary in time due to stochastic fluctuations in environmental variables (Boyce 



et al. 2006). In their simplest guise, models of stochastic demography are the form n t+1 = A f n f 



where n t is the vector of population abundances in the different stages at time t and A t is a 
projection matrix whose entries describe fluxes between stages due to combinations of growth, 
survivorship, and reproduction. These models play a critical role in identifying to what extent 



35 .variation in vital rates alter the stochastic growth rate of a population (Tuljapurkar 1990; Tul 



36 
37 



japurkar et al. 2003 Haridas and Tuljapurkar 2005 Morris et al. 2008). For example, Morris 



et al. (2008) analyzed multiyear demographic data for 36 plant and animal species and found 



that the stochastic growth rate of short-lived species (insects, annual plants, and algae) are 
more negatively affected by increased variation in vital rates than longer-lived species (peren- 
nial plants, birds, ungulates). Correlations between vital rates within a year or autocorrelations 
in vital rates in successive years can significantly effect stochastic growth rates. Analytic ap- 
proximations imply that negative correlations between vital rates can buffer the effects of de- 



mographic stochasticity and thereby increase stochastic growth rates (Tuljapurkar 1990 Boyce 



et al. 2006). In contrast, temporal autocorrelations can increase or decrease stochastic growth 



45 
46 



rates, and increases in these autocorrelations can have larger impacts on growth than increases 



in inter-annual variability of vital rates (Tuljapurkar 1990; Tuljapurkar and Haridas 2006). 
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47 
48 



For populations not experiencing an Allee effect, a positive growth rate at low densities is 
necessary for establishment and range expansion. An important class of models for predicting 
rates of range expansion are integro-difference equation (IDE) models ( Kot et aL]|1996 ) in which 
a dispersal kernel describes the likelihoods that individuals move between locations. For diffusive 
movement, this kernel is a Gaussian and the rate of spread equals 2y/rD where D is the diffusion 



52 rate and r is the intrinsic rate of growth of the population (Kot et al. 1996); the same rate of 



spread first derived by Fisher (1937) for react ion- diffusion equations. While this estimate of 



rate of spread has been applied to many species (Hengeveld 1994), field estimates of dispersal 



55 
56 



kernels typically are leptokurtic and not Gaussian (Kot et al. 1996) and, consequently, this 



earlier work may underestimate invasion speeds. Using stochastic counterparts to these IDE 



57 models, Neubert et al. (2000) showed that serially uncorrected, stochastic fluctuations generate 



58 normally distributed invasion speeds whose variance decays to zero. Hence, invasion speeds may 



59 
60 
61 
62 
63 
64 



jCxhibit unpredictable transients in fluctuating environments. On the other hand, Neubert and 
Caswell ( |2000 ) developed methods to estimate invasion speeds for stage-structured IDE models 



that have been invaluable for identifying how stage-specific vital rates constrain rates of spatial 



spread (Caswell et al. 2003; Jacquemyn et al. 2005 Jongejans et al. 2008). However, the dual 



effects of demographic and temporal heterogeneity on rates of range expansion remains to be 
understood. 



65 Here, we provide a framework for analyzing the simultaneous effects of environmental fluc- 
ee tuations and demographic structure on invasion speeds. We begin by reviewing the work of 



Neubert and Caswell (2000) in constant environments. We extend these models to allow for 



es temporal variation in the projection matrices and dispersal kernels and provide a formula for 

69 asymptotic invasion speeds and normal approximations that describe the variation in these inva- 

70 sion speeds over finite time horizons. Applying these results to two examples, we illustrate how 

71 our results allow one to address questions like "how does the magnitude of variability influence 

72 asymptotic invasion speeds?" and "how do temporal autocorrelations influence the uncertainty 

73 in predicting invasion speeds?" 



Constant environments 



Neubert and Caswell (2000) analyzed invasion speeds for IDE models for stage-structured 



76 populations in constant environments. These models consider structured populations living on 
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77 a continuous one-dimensional habitat and consisting of m stages where n\(x) is the density 

78 of stage i at time t in location x. Let 6 JJ (n 1 (x), . . . , n m (x)) be the contribution of stage j 

79 individuals to stage i individuals at location x. Let k % *(x) be the probability density function 
so for the displacement x moved by an individual transitioning from stage j to i. Under these 
ei assumptions, the dynamics of the population are given by 

/oo m 
k ij (x - y)^{n]{y), . . . , n?{y))4{y) dy 
oo , =1 



82 To simplify the notation, let n t (x) = (n\{x), . . . , n] Ti (x)) / where ' denotes transpose be the vector 

83 of population abundances at time t and location x. Let B(n t (x)) and K(x) denote the m x m 

84 matrices with entries b^(n t (x)) and k^(x), respectively. With this notation, we get the simplified 

85 equation 

/oo 
[K(y - x) oB(n t (y))]n t (y)dy (1) 
-oo 

86 where o denotes the Hadamard product i.e. component wise multiplication. 

87 When the population is unstructured (i.e. m = 1), model ([I]) has traveling waves-solutions 
as that maintain a fixed shape in space and move at a constant speed. If the growth function b(n)n 

89 increases with population density, b(n) is decreasing (i.e. no Allee effect), and the dispersal 

90 kernel possesses a moment-generating function m(s) = J_ k(x)e~ sx dx for < s < s, then the 

91 traveling wave has an asymptotic speed 

c* = min - In \b(0)m(s)] 

0<s<§ s 

92 A population initially concentrated in a finite region of space will never spread faster than 



93 c* and asymptotically will spread at a rate of exactly c* (Weinberger 1982). The linearization 

94 conjecture states that the speed of invasion for a nonlinear model is governed by its linearization 

95 at low population densities, as long as there are no Allee effects and no long-distance density 

96 dependence. This conjecture is supported by many numerical studies. 



97 Relying on the linearization conjecture, Neubert and Caswell (2000) derived a formula for 

98 traveling wave speeds for structured populations. This derivation makes four assumptions: (1) 

99 the matrices B(n) are non-negative and primitive, (2) A = B(0) has a dominant eigenvalue p(A) 

100 that is greater than one, (3) there is negative density dependence i.e. B(n) < An for all n > 

101 where inequalities are taken componentwise, and (4) the kernels (x) have moment generating 
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102 functions m % i(s) defined on some maximal interval < s < s. Under these assumptions, Neubert 

103 and Caswell| ( |2000 ) showed the asymptotic wave speed is given by 



min - logp(A o M(s)) 

0<s<s s 



(2) 



104 where M(s) is the matrix of moment generating functions m l ^{s). 



105 Fluctuating environments 

loe To account for temporal variation in environmental conditions and dispersal rates, we allow 

107 B 4 (n), and K t (x) to depend on time. In which case, 

/oo 
[K t (z - y) o B t (n t (y))} n t (y) dy (3) 
-oo 

108 In order to make use of the linearization conjecture, we place four assumptions on (J3|. First, 

109 the populations exhibit negative density dependence in which case A t := B t (0) < B t (n) for all 
no n > 0. Second, there exists uq > such that A no . . . Ai has all positive entries with probability 
in one. This assumption ensures the effects of the initial distribution of the population disappear in 
112 the limit (demographic weak ergodicity). Third, there exists a s > such that the kernels kf(x) 
in have moment generating functions m % t(s) defined for < s < s. We always choose s to be the 
114 maximum of these values. Without this assumption, solutions may accelerate and have no well- 



H5 defined asymptotic wave speed (Kot et al. 1996). To state the fourth assumption, let M 4 (s) be 



lie the matrix with entries rn t J (s) and H t (s) = M t (s)o A t . We assume that the sequence of matrices, 

in Ho(s), Hi(s), . . . are stationary and ergodic and E[max{ln ||H 4 (s)||, 0}] < oo for all < s < 's. 

us This assumption encompasses many models of environmental fluctuations including periodic, 

n9 quasi-periodic, irreducible Markovian, and auto-regressive models. Under this assumption, the 

120 random version of the Perron- Frobenius theorem ( Arnold et al.||1994 ) implies that for any n > 



121 and w > (i.e. all components are non-negative and at least one is positive) 

lim - ln(Hj_i(s) . . . H (s)n, w) = 7(5) with probability one (4) 

t—>-oo t 

122 where 7(s) is the dominant Lyapunov exponent associated with this random product of matrices 

123 and (•, •) denotes the standard Euclidean inner product. When s = and w = (1, . . . , 1)', this 

124 Lyapunov exponent 

7(0) = lim - ln(A t _i(s) . . . A (s)n, w) 

t— >oo t 
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Fig. 1. — The temporal dynamics of the wave speed t ~ x ° for 250 simulations of the non-linear 
juvenile- adult model. The front of the wave was determined by a threshold of n c = 0.001 with 
equal weight on both stages i.e. w = (1, 1)'. The dashed line is the predicted asymptotic wave 
speed in ([5]). In the inset, a histogram of the waves speeds at t — 500 with the predicted normal 
approximation from the linearization. Parameter values are p = 0, fi = In 40, a = 0.5, a = 1, 
sj = 0.3, and sa = 0.4. 



125 describes the growth rate of the total population size when rare (Tuljapurkar 1990). Since 

126 an invasion can only proceed if the population has the capacity to exhibit growth, our final 
12? assumption is that 7(0) > 0. 

128 When demographic or dispersal rates vary in time, the rate of spatial spread will vary in 

129 time. To quantify this rate of spread, assume there is a weighting w > of the different stages 

130 and a weighted population abundance, n c > 0, above which the population is observable. At 

131 time t, let X t denote the location furthest from the invasion's origin at which (n^s^w) > n c . 
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T=50 T=100 




Fig. 2. — Spatio-temporal dynamics of range expansion for the juvenile-adult model. Spatial 
distribution and abundance of juveniles (in shaded red) and adults (in shaded grey) plotted at 
the indicated times. Parameters (p = 0, \i = In 40, a = 0.1, a — 1, sj — 0.3, and sa = 0.4) are 
such that local dynamics are chaotic. 



132 The average speed Ct of the invasion by time t is given by Xt ~ x ° where xq is the introduction 

133 site. In Appendix A we show that with probability one, Xt ~ x ° is asymptotically bounded above 

134 by 

c = mm (5) 

0<s<S s 

135 whenever rio(x) has compact support. The linearization conjecture and our numerical results 

136 (Fig. [T]) suggest that c* is not only an upper bound but in fact equals the asymptotic wave speed 

137 with probability one. In the special case of periodic environmental fluctuations (i.e. there exists 

138 a natural number p such that A t + P = A t and K f+P (x) = K t (x) for all t and x), the asymptotic 

139 invasion speed is given by 



c* = min — lnp ( TT A t M t (s) ) 
s>0 sp ^±± J 



wo where p denotes the dominant eigenvalue of a matrix. 



i4i n Under the additional assumption that matrices H t (s) are rapidly mixing (see, e.g., |Heyde 



and Cohen 



(1985)), Appendix B shows that there exists a > such that the average speeds Xt ~ x ° 

143 are asymptotically normal with mean c* and standard deviation a/(s*\/i). This approximation 

144 is consistent with numerical simulations of the nonlinear models (inset of Fig. [I]). In the special 

145 case of an unstructured population m — 1, this result extends Neubert et al. (2000)'s work on 
we uncorrelated environments to correlated environments. In this case, the central limit theorem 
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for stationary sequences of random variables (Durrett 1996) implies 



a 2 = Var[lnHi(s*)] + 2 ^ Cov[ln Hi(s*), In H 4 (s*)] 



(6) 



i=2 



where s* is such that c* = j(s*)/s*. For unstructured populations, equation ^ implies that pos- 
itive temporal autocorrelations increase the variability in asymptotic wave speeds while negative 
autocorrelations reduce this variability. In contrast, as 7(5) = E[ln ifi(s)], temporal autocorre- 
lations have no effect on the mean c* invasion speed for unstructured populations. 

To illustrate the applicability of our results to structured populations, we consider two 
examples. The first example is a simple juvenile-adult structure model with continuous variation 
in fecundity rates. We use this example to illustrate the accuracy of the linear approximation, the 
dynamics of spatial spread, and the effects of noise amplitude and color on invasion speeds and 
stochastic growth. In the second example, we analyze the effects of environmental stochasticity 
on invasion speeds for an herbaceous perennial herb, Calathea ovandensis (Marantaceae) . This 



158 example was studied for constant environments by Neubert and Caswell (2000). Here, we 



159 
160 



illustrate the effects of the frequency of poor environmental conditions and temporal correlations 
on invasion speeds and their uncertainty using empirical data. 



161 Example: Fluctuating fecundities. We consider a population with two stages, juvenile 

162 and adult, in which dispersal occurs following reproduction. We assume that the per-capita 

163 fecundity f t exp(—an 2 ) is density-dependent where n 2 is the adult density at time t, f t is the 

164 maximal fecundity at time t, and a > measures the intensity of intraspecific competition. Let 

165 sj denote the fraction of juveniles surviving to adulthood, and sa denote the fraction of adults 
lee that survive to the next year. Under these assumptions, 

/0 f t eM-an 2 )\ 

\SJ S A J 

167 We assume that log/* are normally distributed with mean /i, variance a 2 , and temporal corre- 

168 lation p between log/t and log/t+i. Dispersal follows birth with a Laplacian distribution with 

169 variance 2b 2 . Hence, the moment generating matrix for this model is 



M(s) = 
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Fig. 3. — The stochastic growth rate 7 and invasion speed c* plotted as a function of the 
standard deviation of the log fecundity log/ t . Length of error bars correspond to the standard 
deviation in the growth rate after 400 generations. The blue, white, and red curves correspond 
to temporal correlations of p = —0.5,0 and 0.5. Parameters are such that E[/J = 4, sj = 0.3, 
and sa = 0.4. Simulations for estimates ran for 10,000 generations. 

170 which is defined for < s < s = 1/6. When adults are highly fecund, the interaction between 

171 locally unstable population dynamics and environmental stochasticity can generate complex 

172 spatial-temporal patterns of range expansion (Fig. [2]). None the less, simulations of the full 

173 nonlinear mode suggest that linearization conjecture holds (Fig. [TJ. 

174 To understand the effects of temporal variation and correlations on wave speeds and popula- 

175 tion growth, we computed 7(0), c*, and standard deviations of both of these quantities for a range 

176 of /i and a 2 values. In all of these simulations, we held the mean fecundity E[/t] = exp(yU — <r 2 /2) 

177 constant at 4. The simulations show that the population growth rate 7(0) and the asymp- 

178 totic wave speed decreases both with increasing temporal variation and increasing temporal 

179 autocorrelations (Fig. [3]). In particular, when this variation is too high, the population has a 

180 negative stochastic growth rate and does not propagate. Autocorrelations have little effect on 

181 the variability of the stochastic growth rates, but variation in invasion speeds increases with 

182 autocorrelations. 
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183 Example: Calathea ovandensis. Calathea ovandensis is an understory monocot, found 



184 in neotropical lowland rain forests (Horvitz and Schemske 1995). Deciduous during the dry 



185 
186 
187 
188 



season, C. ovandensis plants reinitiate growth from rhizomes during the rainy season and bear 



fruit capsules, which dehisce to expose large seeds that fall to the forest floor (Horvitz and 
Beattie 1980[ ). Seeds are myrmecochorous (adapted for ant dispersal), bearing oily white arils 



(eliaosomes) that are considered arthropod prey mimics (Carroll and Janzen 1973). Ponerine 
ants carry seeds to their nests like prey, remove and consume the arils, and bury seeds in 



190 nitrogen-rich garbage piles of dead insect prey near the nest entrance (Horvitz and Beattie 



1980 Horvitz and Schemske 1995). 



193 
194 



At Horvitz and Schemske's study site in southern Mexico, four ant species act as disperal 
agents: Pachycondyla apicalis, P. harpax, Solenopsis geminata, and Wasmannia auropuncata 
( Horvitz and Schemske|19"95 ). Calathea ovandensis do not propagate vegetatively, and empirical 



195 distributions of seedlings are well matched to ant dispersal distances (Horvitz and Schemske 



196 
197 



1986). Horvitz and Schemske (1995) measured demographic rates, constructing 16 projection 



matrices for four sites over four years, and found that C. ovandensis' highest population growth 



198 rates occurred during el Nino years and in plots affected by treefall gaps (Horvitz and Schemske 



1995 C. Horvitz and Caswell 1997). 



We selected 3 stage-structured matrices that represent 1) the highest population growth, 



2oi observed during an el Nino year (A = 1.2477; plot 2, 1982-83; also used in Neubert and Caswell 



202 
203 



(2000)), 2) lowest (negative) growth in the same plot (A = 0.9051; plot 2; 1984-85), and 3) the 



worst year-plot combination overall (A = 0.7356; plot 3, 1984-85). Individuals are classified as 
seeds, seedlings, juveniles, pre-reproductives, and small, medium, large and extra-large adult 



205 plants (for life cycle graph, see Neubert and Caswell (2000)). Seed dispersal involves transitions 



206 from most stages. We follow Neubert and Caswell (2000) in assuming a Laplace distribution 



of dispersal distances of four ant species and weighting dispersal by the relative abundance of 
each species ( Horvitz and Schemske||T995 Neubert and Caswell||2000 ). This treatment assumes 



that seed fates are independent of dispersal agent (Neubert and Caswell 2000), and yields the 
dispersal kernel matrix 



K(x -y) = ^PiKi 



x 



y) 



i=l 



where pi are the proportion of seeds dispersed by each ant species, and Ki(x — y) is a matrix 
whose first row contains Laplace dispersal kernels with generating functions m,i(s) = 1/(1 — b 2 s 2 ) 
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213 with parameter b{, representing the mean distance a seed is moved by a particular ant species 



2u as reported by Horvitz and Schemske (1986), and Delta functions for the remaining entries 



215 (Neubert and Caswell 2000). The moment-generating function matrix is 

4 

i=i 

216 where Mj(s) is the moment generating function for K^x). 

217 We ran two types of simulations representing (1) random fluctuations between environmen- 

218 tal states corresponding to the highest population growth and the lowest population growth in 

219 the same plot (i.e. randomly varying between the matrices with dominant eigenvalues A = 1.2477 

220 and A = 0.9051) and (2) random fluctuations between environmental states corresponding to the 

221 highest and lowest overall population growth across all sites (i.e. A = 1.2477 and A = 0.7356). 

222 For each of these simulations, switches between "good" and "bad" years were modeled by a 

223 Markov chain in which the probability of transition from a bad year to a good year is q and the 

224 probability of transition from a good year to a bad year is p. For this model of the environmen- 

225 tal dynamics, the long-term frequency of good years is and the temporal autocorrelation 

226 between environmental states is 1 — q — p. For both scenarios, simulations were run for all 

227 combinations of p and q values to determine the combined effects of frequency of good years 

228 and temporal correlations on stochastic growth and invasion speed. 

229 Both the stochastic growth rate and invasion speed were strongly correlated with the overall 

230 frequency of good years, while correlations between environmental states had a relatively small 

231 effect on asymptotic growth rates and invasion speeds (Fig. [4]). However, positive correlations 

232 between environmental states did increase growth rates and invasion speeds, while negative 

233 correlations reduced both (left and middle column of Fig. [1]). For growth rates, the effect of 

234 correlations was greatest when the frequency of good and bad years were equal {q/{q + p) = 0.5) 

235 (left column of Fig. [4]), presumably since the largest number of combinations of (p, q) values 

236 can generate this frequency. Environmental correlations had the greatest impact on invasion 

237 speeds near the threshold between positive and zero growth (middle column of Fig. [4]), where 

238 invasion speeds also exhibited the greatest variance (right column of Fig. [4]). As a result, near the 

239 threshold in parameter values between positive and negative growth rates (and hence positive 

240 or zero invasion speed), the magnitude and sign of correlations between environmental states 

241 determined whether populations grew and spread or declined. For example, for scenario one 
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Top row: ^=1.2477, X 2 =0.9051 
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Bottom row: ^=1.2477, X 3 =0.7356 
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Fig. 4. — Stochastic growth rates and invasion speeds (c*) for Calathea ovandensis. Stochastic 
growth rates (left column), asymptotic invasion speeds (middle column), and standard deviation 
in invasion speeds after 400 years (right column) are plotted as a function of the frequency q/(q + 
p) of good years for all probabilities p and q of switching from "good" to "bad" years and vice 
versa. Temporal correlations associated with the different p, q combinations are plotted in color 
as indicated in the legend. Figures in the top row correspond to random fluctuations between 
environmental states corresponding to the highest population growth and the lowest population 
growth in the same plot (i.e. A = 1.2477 and A = 0.9051). Figures in the bottom row correspond 
to random fluctuations between environmental states corresponding to the highest and lowest 
overall population growth across all sites (i.e. A = 1.2477 and A = 0.7356). Simulations ran for 
100,000 generations. 

242 (top row of Fig. [4]), strong positive correlations allowed the population to grow when only 35% 

243 of years were good, while strong negative correlations generated negative growth rates under 

244 the same scenario even at nearly even proportions of good and bad years. 
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245 Like growth rates and invasion speed, variance in invasion speed was also strongly affected 

246 by both environmental correlation and the frequency of good years (right column of Fig. [3]) . Most 

247 notably, invasion speed variance increased approaching the threshold of zero population growth 

248 rates in both scenarios, though most dramatically in the second (bottom row of Fig. [1]), and 

249 declined as the proportion of good years increased. Variance in invasion speed also increased 

250 as environmental correlations shifted from strongly negative to strongly positive. Hence, the 

251 combinations of environmental states that yielded the fastest invasion speeds also showed the 

252 highest variance in invasion speed. 

253 Comparing between scenarios, the scenario with a higher average growth rate between 

254 years (top row of Fig. [3]) exhibited a lower threshold for positive growth rates with respect to 

255 the frequency of good years and overall lower variance in invasion speeds, as expected. Since 

256 positive invasion speeds increased roughly linearly with the proportion of good years, a necessary 

257 outcome was that each increment of change in the proportion of good years led to a greater 

258 proportional increase in invasion speed in the "worse" scenario (bottom row of Fig. [4]) than 

259 the scenario in which average growth rate across years was higher (top row of Fig. [4]). In the 

260 context of real systems, this means that we might expect to observe particularly steep increases 

261 in invasion speed in response to increased frequency of good years in regions where bad years 

262 are particularly bad. In particular, changes in the frequency of El Nino years or in the temporal 

263 correlation between El Nino events, could substantially affect rates of spread of Calathea. 



Discussion 



How environmental variation influences population dynamics and biodiversity is a central 



266 question in ecology (Andrewartha and Birch 1954 May 1975 Tuljapurkar 1982 Chesson 2000 



Holyoak et al. 2005). In our era of increasing climate variability and accelerating global trans- 



268 port of species at rates and distances far beyond their innate capacity (Vitousek et al. 1997), 



269 this question is not only academic but also directly relevant to wildlife management and the 

270 maintenance of regional biodiversity. The ability to forecast and mitigate impacts of success- 

271 ful invaders, and changes in the distribution of native species, relies on an enhanced under- 

272 standing of the potential impacts of environmental variability on population growth rates and 



273 invasion speeds (Neubert et al. 2000; Collingham and Huntley 2000). By extending a class of 



stage-structured integrodifference equation models to variable environments, we developed and 
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275 applied a methodology for examining how environmental variability and temporal correlations 

276 impact stochastic growth rates, invasion speeds, and uncertainty in invasion speed predictions. 

277 For structured populations in fluctuating environments, our analysis reveals that invasion 

278 speeds are approximately normally distributed and identifies data requirements for estimating 



279 the mean and variance of invasion speeds. As shown by Neubert et al. (2000) for unstructured 

280 populations in serially uncorrelated environments, we found that variance in invasion speeds de- 

281 cay in time. Hence, the greatest uncertainty in predicting the spatial extent of an invader occurs 

282 in the earliest stages of spatial spread. The data requirements for quantifying this uncertainty 

283 are two- fold. To capture the effects of environmental fluctuations on demography, estimates of 

284 survivorship, fecundity, and transition rates between different stages are needed for the range of 

285 environmental conditions most likely to be encountered by the population ( Caswell||2001 Boyce 



et al. 2006). Our methods also require estimates for the dispersal kernel associated with each 



287 demographic transition under an appropriate range of environmental conditions. While getting 

288 such estimates is challenging, our methods can help focus empirical efforts by identifying which 

289 stage-specific demographic rates have the largest effect on the mean and variance of the invasion 

290 speed. More generally, they can be used as an exploratory tool to identify key features of how 

291 environmental variability and stage-structure interact to determine invasion speeds. 

292 Our examples illustrate several important general effects of environmental variability on 

293 invasion speeds, and point to the ways in which the outcomes of environmental-demographic 

294 interactions might vary between species. The results of our two-stage example with fluctuating 

295 fecundities showed growth rates and invasion speeds declining as the magnitude of variance 

296 in vital rates increased. This result is consistent with earlier work showing that increases in 

297 the magnitude of environmental (e.g. climate) variability tend to reduce population growth 

298 rates ( Lewontin and Cohen|1969 ). Theory and empirical findings predict that the magnitude of 

299 these effects on population growth will depend on a species' life history: short-lived species are 



more negatively affected by increased variation in vital rates than longer-lived species ( Chesson 



1985 Morris et al.||2008 ). Although not investigated here, it seems reasonable to conjecture that 



302 these life history traits will have similar effects on invasion speeds. 

303 While increasing environmental variability reduces growth and invasion speed, correlations 

304 between environmental states can have a substantial impact on the direction and variance in 

305 growth and wave speed. Our analysis of a general unstructured model and two structured models 
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shows that positive temporal correlations between environmental states (reddened spectrum in 
figures) created greater variability in growth rates and invasion speeds, even at low levels of envi- 
ronmental noise. This finding is consistent previous work, which shows potentially strong effects 



310 
311 



pf positive temporal correlations in vital rates on variance in population growth (Tuljapurkar 
and Orzack||1980 Tuljapurkar ||1982" Runge and Moen"||1998 ). In contrast to the consistent effect 



of correlations on variance, the direction of effects on expected growth rates and invasion speeds 
differed between examples. The general unstructured model shows no effect, which suggests 
that demographic structure allows for effects to go in either direction. Positive environmental 
correlations yielded negative effects on both stochastic growth rates and invasion speeds in the 
two-stage model, but positive effects in the empirically based Calathea model. Similarly, the 
direction of effects of correlations have been shown to vary in a life history context as well: serial 
correlations can have strongly positive or negative effects on fitness, depending on the structure 
of life history and the nature of within-year (among reproductive classes) and between-year 



319 correlations in fertilities in response to environmental variation (Tuljapurkar et al. 2009). 



We also found that the magnitude of effect of environmental correlations differed between 
examples. In the Calathea example, the effect of correlations was small relative to the impor- 
tance of the overall proportion of years with positive growth rates, and was more limited than 
in the two-stage model. Limited effects of temporal environmental correlations on long-term 



324 stochastic growth rates have been observed in other empirical studies (Silva et al. 1991), and 



325 theoretical studies have found relatively minor effects of temporal correlation in vital rates on 



326 stochastic growth rates (Tuljapurkar 1982 Fieberg and Ellner 2001). However, the effects of 



327 temporal autocorrelation can be great under some circumstances, depending on life history and 



how well vital rates "remember" the past (Tuljapurkar and Haridas 2006). Populations with 



weak demographic damping (strong transient effects) tend to accumulate greater variability 
over time than those with strong damping, and these effects vary with the magnitude of serial 



correlation (Tuljapurkar and Haridas 2006). 



333 
334 



Although the Calathea example showed relatively small effects of correlations overall, cor- 
relations still had the potential to substantially influence growth rates and invasion speed under 
some conditions, especially for populations near the threshold of positive and negative growth 
rates. In these threshold cases, strongly positive correlations resulted in the fastest invasion 
speeds, and also generally the greatest variance in those speeds. This suggests that, on the 
ground, changes in environmental correlations could lead to sudden or episodic changes in inva- 
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338 sion rates for species experiencing threshold growth conditions. These kinds of dynamics might 

339 contribute to common observations of a lag between establishment and spread of most introduced 

340 species, even those that eventually become invasive, alongside other factors such as adaptation 

341 |and hybridization with native species ( |Ellstrand and Schierenbeck||2000| |Mack et al.|2000HSakai 



et al. 2001). Among native species responding to climate change, these results suggest that 



343 we might expect to observe great variation in rates of range expansion, since populations at 

344 the edge of species ranges are most likely to exist close to the threshold between positive and 

345 negative population growth. In either case, the increased variability in stochastic growth rates 

346 and invasion speed around threshold conditions will make biological expansions and invasions 

347 harder to forecast. For species characterized by substantial interannual variation and volatility 

348 in vital rates, accurately forecasting rates of spread may be particularly challenging. 

349 Several factors must be taken into careful account in applying integrodifference equation 

350 models to invasions on real landscapes. For one, it is critical to understand the basic nature 

351 of relevant environmental factors and whether a given form of disturbance imposes different 



352 consequences when rare than it does when chronic (Tuljapurkar and Haridas 2006). For example, 

353 Tuljapurkar and Haridas point out that an increase in the frequency of hurricanes may weaken 

354 autocorrelation by preventing successional change, thereby effectively reducing demographic 

355 variability and the effects of autocorrelation over time despite the outward perception of dramatic 

356 environmental variation. Similarly, single years or short-term strings of conditions (e.g., warm 

357 years) may boost demographic rates for a given species, while extended autocorrelated series 

358 would lead to more substantial ecological changes (e.g., aridification) that negatively impact 

359 habitat suitability and represent non-Markovian dynamics. The models we present are applicable 

360 to increasing climatic variability but not to sustained directional changes over time. In the 

361 context of climate change, these models are likely sufficient to describe many regions of the globe, 

362 where increasing climate variability is expected to obscure the effects of concurrent directional 

363 trends over much of the next century (IPCC 2007), but not appropriate for regions where 

364 trending changes in temperature are more acute. 

365 Despite these limitations, our methods provide a first step in estimating invasion speeds 

366 for populations experiencing stage-specific temporal variation in survivorship, transition rates 

367 between stages, fecundity, and dispersal. Using these methods to investigate the effects of 

368 temporal correlations on rates of invasion across a wider range of life histories and dispersal 

369 types will be an important next step. 



-17- 



syo REFERENCES 

371 H.G. Andrewartha and L.C. Birch. The distribution and abundance of animals. University of 

372 Chicago Press, 1954. 

373 L. Arnold, V. M. Gundlach, and L. Demetrius. Evolutionary formalism for products of positive 

374 random matrices. Annals of Applied Probability, 4:859-901, 1994. 

375 M. S. Boyce, C. V. Haridas, C. T. Lee, and the NCEAS Stochastic Demography Working Group. 

376 Demography in an increasingly variable world. Trends in Ecology & Evolution, 21 (3): 141 

377 - 148, 2006. 

378 D. Schemske C. Horvitz and H. Caswell. The importance of life history stages to population 

379 growth: prospective and retrospective analyses. In S. Tuljapurkar and H. Caswell, editors, 

380 Structured population models in marine, terrestrial and freshwater systems, pages 247- 

381 272, New York, New York, USA, 1997. Chapman and Hall. 

382 H. Caswell. Matrix Population Models. Sinauer, Sunderland, Massachuesetts, 2001. 

383 H. Caswell, R. Lensink, and M. G. Neubert. Demography and dispersal: Life table response 

384 experiments for invasion speed. Ecology, 84(8): 1968-1978, 2003. 

385 P.L. Chesson. General theory of competitive coexistence in spatially-varying environments. 

386 Theoretical Population Biology, 58:211-237, 2000. 

387 P.L. Chesson. Coexistence of competitors in spatially and temporally varying environments: 

388 a look at the combined effects of different sorts of variability. Theoretical Population 

389 Biology, 28:263-287, 1985. 

390 Y.C. Collingham and B. Huntley. Impacts of habitat fragmentation and patch size upon migra- 

391 tion rates. Ecological Applications, 10(1):131-144, 2000. 

392 R. Durrett. Probability: Theory and examples. Duxbury Press, Belmont, CA, 1996. 

393 D.R. Easterling, G.A. Meehl, C. Parmesan, S.A. Changnon, T.R. Karl, and L.O. Mearns. Cli- 

394 mate extremes: observations, modeling, and impacts. Science, 289:2068, 2000. 



- 18 - 



395 N. C. Ellstrand and K. A. Schierenbeck. Hybridization as a stimulus for the evolution of inva- 

396 siveness in plants? Proceedings of the National Academy of Sciences of the United States 

397 of America, 97:7043-7050, 2000. 

398 C. S. Elton. The Ecology of Invasion by Animals and Plants. Methuen, 1958. 

399 J. Fieberg and S.P. Ellner. Stochastic matrix models for conservation and management: a 
wo comparative review of methods. Ecology Letters, 4(3):244-266, 2001. 

401 R.A. Fisher. The wave of advance of advantageous genes. Ann. Eugenics, 7:269, 1937. 

402 C. V. Haridas and S. Tuljapurkar. Elasticities in variable environments: properties and impli- 

403 cations. The American Naturalist, 166:481-495, 2005. 

404 R. Hengeveld. Small-step invasion research. Trends in Ecology & Evolution, 9:339-342, 1994. 

405 C.C. Heyde and J.E. Cohen. Confidence intervals for demographic projections based on products 

406 of random matrices. Theoretical Population Biology, 27(2):120-153, 1985. 

40? M. Holyoak, M. A. Leibold, and R. D. Holt, editors. Metacommunities: Spatial dynamics and 

408 ecological communities. University of Chicago Press, 2005. 

409 C. C. Horvitz and A. J. Beattie. Ant dispersal of Calathea (Marantaceae) seeds by carnivorous 

410 ponerines (Formicidae) in a tropical rain forest. American Journal of Botany, 67:321-326, 

411 1980. 

412 C. C. Horvitz and D. W. Schemske. Seed dispersal of a neotropical myrmecochore: variation in 

413 removal rates and dispersal distance. Biotropica, 18:319-323, 1986. 

414 C. C. Horvitz and D. W. Schemske. Spatiotemporal variation in demographic transitions of a 

415 tropical understory herb: projection matrix analysis. Ecological Monographs, 65:155-192, 
1995. 

417 H. Jacquemyn, R. Brys, and M. G. Neubert. Fire increases invasive spread of molinia caerulea 

418 mainly through changes in demographic parameters. Ecological Applications, 15(6):2097- 
«9 2108, 2005. 



- 19 - 



420 E. Jongejans, K. Shea, O. Skarpaas, D. Kelly, A.W. Sheppard, and T.L. Woodburn. Dispersal 

421 and demography contributions to population spread of Carduus nutans in its native and 

422 invaded ranges. Journal of Ecology, 96:687-697, 2008. 

423 C.S. Kolar and D.M. Lodge. Progress in invasion biology: predicting invaders. Trends in Ecology 

& Evolution, 16(4):199-204, 2001. 

425 M. Kot, M. A. Lewis, and P. van den Driessche. Dispersal data and the spread of invading 

426 organisms. Ecology, 77:2027-2042, 1996. 

427 R. C. Lewontin and D. Cohen. On population growth in a randomly varying environment. 



428 Proceedings of the National Academy of Sciences of the USA, 62:1056-1060, 1969. 

429 R.N. Mack, D. Simberloff, W. Mark Lonsdale, H. Evans, M. Clout, and F.A. Bazzaz. Biotic 

430 invasions: causes, epidemiology, global consequences, and control. Ecological applications, 

431 10(3):689-710, 2000. 

432 R. M. May. Stability and Complexity in Model Ecosystems, 2nd edn. Princeton University Press, 

433 Princeton, 1975. 

434 W.F. Morris, C.A. Pfister, S. Tuljapurkar, C.V. Haridas, C.L. Boggs, M.S. Boyce, Emilio M. 

435 Bruna, Don R. Church, Tim Coulson, D.F. Doak, S. Forsyth, J. Gaillard, C.C. Horvitz, 

436 S. Kalisz, B.E. Kendall, T.M. Knight, C.T. Lee, and E.S. Menges. Longevity can buffer 

437 plant and animal populations against changing climatic variability. Ecology, 89:19-25, 

438 2008. 

439 M. G. Neubert and H. Caswell. Demography and dispersal: Calculation and sensitivity analysis 

440 of invasion speed for structured populations. Ecology, 81:1613-1628, 2000. 

441 M. G. Neubert, M. Kot, and M. A. Lewis. Invasion speeds in fluctuating environments. Pro- 

442 ceedings of Royal Society of London B, 267:1603-1610, 2000. 

443 M.C. Runge and A.N. Moen. A modified model for projecting age-structured populations in 

444 random environments. Mathematical Biosciences, 150(1) :21-41, 1998. 

445 A. K. Sakai, F. W. Allendorf, J. S. Holt, D. M. Lodge, J. Molofsky, K. A. With, S. Baughman, 

446 R. J. Cabin, J. E. Cohen, N. C. Ellstrand, D. E. McCauley, P. O'Neil, I. M. Parker, J. N. 



- 20 - 



44? Thompson, and S. G. Weller. The population biology of invasive species. Annual Review 

448 of Ecology and Systematics, 32:305-332, 2001. 

449 J. F. Silva, J. Raventos, H. Caswell, and M. C. Trevisan. Population responses to fire in a 

450 tropical savanna grass, andropogon semiberbis: a matrix model approach. Journal of 



451 Ecology, 79:345-355, 1991. 

452 C. Tebaldi, K. Hayhoe, J.M. Arblaster, and G.A. Meehl. Going to the extremes. Climatic 

453 Change, 79(3):185-211, 2006. 

454 S. Tuljapurkar. Population dynamics in variable environments. II. Correlated environments, 

455 sensitivity analysis and dynamics. Theoretical Population Biology, 21:114-140, 1982. 

456 S. Tuljapurkar. Population Dynamics in Variable Environments. Springer- Verlag, New York, 

1990. 

458 S. Tuljapurkar and C. V. Haridas. Temporal autocorrelation and stochastic population growth. 

459 Ecology Letters, 9:327-337, 2006. 

460 S. Tuljapurkar, C.C. Horvitz, and J.B. Pascarella. The many growth rates and elasticities of 

461 populations in random environments. The American Naturalist, 162(4):489-502, 2003. 

462 S. Tuljapurkar, J.M. Gaillard, and T. Coulson. From stochastic environments to life histories 

463 and back. Philosophical Transactions of the Royal Society Part B, 364:1499-1509, 2009. 

464 S.D. Tuljapurkar and S.H. Orzack. Population dynamics in variable environments I. Long-run 

465 growth rates and extinction. Theoretical Population Biology, 18(3):314-342, 1980. 

466 P.M. Vitousek, H.A. Mooney, J. Lubchenco, and J.M. Melillo. Human domination of Earths 

467 ecosystems. Science, 277:494-499, 1997. 

468 H.F. Weinberger. Long-time behavior of a class of biological models. SI AM Journal on Mathe- 

469 matical Analysis, 13:353, 1982. 



This preprint was prepared with the AAS IAT^X macros v5.2. 



21 



470 A. Mathematical Details 

471 We begin by considering solutions to the linear equation 

/oo 
[K t (x-y) o A t ] n t (y) dy 
-oo 

472 Consider no(x) = ue~ sx with u > 0. We claim that 

n t (x) =11^(3)... H (s)ue- sx 

473 Indeed, we show this inductively. For t — 1, 

/oo 
[K (y-x)oA ] ne-^dy 
-oo 

/oo 
[K (z)oA ] e- sz dzue~ sx 
-oo 

: H ( S )ue— 



ndxl 



where z = y — x 



474 Next we assume that (A. 2) holds for £ — 1. Then 



/oo 
[K t (?/ - x) o A t ] n t (y) 
-oo 

/oo 
[K t (y - x) o A,] H^s) . . . H ( S )ue-^ dy 
-oo 



oo 
oo 



(Al) 



(A.2) 



[K t (z) o A t ] e~ sz dz H t _ x (a) . . . H (s)u e - 
H t ( S )H t _ 1 ( S )...H (s)ue- s * 



where z = y — x 



475 which completes the inductive proof of (A.2). 



476 Let X t (s) be such that (n t (X t (s)), w) = n c . Then (n (X (s)), w) = (n t (X t (s)), w) and 

(ue- sX °W,w) = (H t _ 1 ( S )...H ( S )ue- sX *W,w) 



477 (\A.2y imply that 



478 equivalent ly 



Y f \ y ( % 1 , (H f _i(s)...H (s)u,w) 

A t (s) - A (s) = - In , r 

s (u, w) 



479 Equation (HI) implies that lim 



t— »oo 



X t (s)-X ( S ) _ 7 (s) 



with probability one. 
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480 Now consider an initial condition n (x) with compact support for the non-linear model ([3]). 

481 Given any s > 0, choose u > such that n (x) < ue~ sx . Our assumption that B t (n) < A f 

482 implies that 

n t (x) <H t _ 1 (s)...H (s)ue- sa: 

483 Thus X t < X t (s) where X t is such that (n t (X t ), w) = n c . Hence, with probability one, 

limsup — - — — < ry(s)/s 

t— >oo t 

484 Since this holds for any < s < s, it follows that 

limsup — ^< min j(s)/s 

t^oo t 0<s<s 

485 with probability one. The Linearization Conjecture implies that 

X t — Xo 

hm 



min r y(s)/s =: c* 

t^oo t 0<s<§ 



486 B. Log-normal approximation 

48? Let s* be such that r y(s*)/s* = c*. Assume that Hj.(s*), H 2 (s*), . . . are rapidly mixing (see, 



488 ^.g., Tuljapurkar (1990); Heyde and Cohen (1985) for definitions). Then Theorem 1 in (Heyde 



and Cohen 1985 ) implies that there exists a > such that 

ln(H t _i(s*)...H (s*)u, w > ~c*s*t 
Via 



490 converges in distribution to a standard normal as t — > oo. 



